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ABSTRACT 


Several multigrid schemes are considered for the numerical computation of viscous 
hypersonic flows. For each scheme, the basic solution algorithm employs upwind spatial 
discretization with explicit multistage time stepping. Two-level versions of the various 
multigrid algorithms are applied to the two-dimensional advection equation, and Fourier 
analysis is used to determine their damping properties. The capabilities of the multi- 
grid methods are assessed by solving two different hypersonic flow problems. Some new 
multigrid schemes, based on semicoarsening strategies, are shown to be quite effective in re- 
lieving the stiffness caused by the high-aspect-ratio cells required to resolve high Reynolds 
number flows. These schemes exhibit good convergence rates for Reynolds numbers up to 
200 x 10 6 . 


INTRODUCTION 


In the past several years, multigrid has been used to accelerate the convergence of 
Navier-Stokes computations for a variety of flow problems at both subsonic and transonic 
speeds (refs. 1 and 2). More recently, multigrid methods with either central or upwind 
differencing have been applied to viscous hypersonic flows to achieve convergence rates 
that approach those obtained at lower Mach numbers and moderate Reynolds numbers 
(Re < 10 7 ). However, at the higher Re values experienced by high-speed flight vehicles, 
a dramatic slowdown occurs in the convergence rate. One reason for this slowdown is the 
deterioration in the high-frequency damping of the multigrid driving scheme caused by the 
very high-aspect-ratio cells that occur in the computational mesh in order to resolve the 
thin boundary layers. 

The present paper describes an effort to understand and improve the use of multigrid 
schemes for the computation of viscous hypersonic flows. First, various two-level multigrid 
schemes both with and without semicoarsening axe introduced. Then we use a Fourier 
analysis of the schemes, applied to the two-dimensional convection equation, to reveal the 
behavior of their components. For each multigrid approach, the solver uses an upwind dis- 
cretization combined with an explicit multistage scheme. We next consider the numerical 
solution of the Navier-Stokes equations for hypersonic flows. The basic elements of the 
flow solver for these equations are summarized. Some details concerning the application of 
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the time-stepping scheme to fine- and coarse-grid problems are presented. The extension 
of the two-level schemes to multilevel ones is then discussed. Elements of multigrid that 
are of particular importance for high-speed flow computations are given. In the results 
section, we consider two different hypersonic flow problems to assess the capabilities of the 
multigrid schemes. Semicoarsening is shown to be quite effective in relaxing the stiffness 
that arises from the resolution of thin boundary layers. 

MULTIGRID METHOD AND STRATEGIES 

The multigrid approach is based on the full approximation scheme of Brandt (ref. 3). 
The grid transfer operators are those considered by Jameson (ref. 4). Coarser meshes 
are obtained by eliminating alternate mesh points in each coordinate direction. Both the 
solution and the residuals are restricted from fine to coarse meshes. A forcing function 
is constructed so that the solution on a coarse mesh is driven by the residuals collected 
on the next finer mesh. The corrections obtained on the coarse mesh are interpolated 
back to the fine mesh. The multigrid schemes investigated within the present work are 
displayed in Figure 1. Figure 1(a) shows a two-level scheme with full coarsening. Re- 
striction of the solution from the fine mesh (m,n) to the coarse mesh (m/2,n/2) is done 
by injection, whereas full weighting is used for the restriction of the residuals. Prolonga- 
tion of the corrections is done by bilinear interpolation. Figure 1 (b) shows a scheme with 
semicoarsening in the different coordinate directions. Again, injection and full weighting 
are used in the restriction process. The corrections obtained on the coarse meshes are 
averaged before they are added to the current fine mesh solution which is indicated by the 
numbers at the “up” arrows. Because of this averaging, half of the individual corrections 
on the coarse meshes are lost. We, therefore, anticipate that the scheme in Figure 1(a) 
should be computationally more efficient, provided that enough high-frequency damping 
ran be obtained with the smoothing scheme of the fine mesh. In order to overcome this 
deficiency of the semicoarsening scheme, two more variants are considered. For the scheme 
of Figure 1(c), the solutions on the coarse meshes are computed sequentially. Hence, the 
corrections obtained on the (m/2,n) mesh can be used to update the (m,n/2) mesh before 
time stepping (as indicated by the horizontal arrow). The sequential update of the second 
coarse mesh allows the full corrections to be passed up to the fine mesh. Note that this 
multigrid variant is not compatible with the idea of parallel computations. An interesting 
compromise between the schemes of Figures 1(b) and 1(c) was suggested by Van Rosendale 
based on the work of ref. 5 (Figure 1(d)). Here, only the corrections common to both 
of the coarse meshes, (m/2,n) and (m,n/2), are averaged, whereas the corrections to the 
modes that live either on (m/2,n) or on (m,n/2) are passed to the fine mesh in full. This 
scheme does allow parallel computations for the coarse meshes. 

FOURIER ANALYSIS OF THE SCALAR ADVECTION EQUATION 

A crucial factor in constructing an effective multigrid method is the selection of a 
smoothing or driving scheme. Local mode (Fourier) analysis is generally applied to evaluate 
possible smoothers on the basis of stability and high-frequency damping properties. The 
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screening of schemes is often performed with a single-grid analysis. Since a stable single- 
grid scheme may not be stable for the multigrid process, the behavior of a smoother with 
a particular multigrid strategy is needed. In addition, the multigrid process can have 
a substantial impact on the performance of the multigrid method. In fact, as we will 
demonstrate in this paper, semicoarsening can provide significant improvement, relative to 
full coarsening, in the damping of the multigrid, especially when a strong mesh anisotropy 
is present due to the high-aspect-ratio cells. 

In ref. 4, Jameson models a multigrid scheme as a multilevel uniform scheme and 
analyzes the stability of this scheme when applied to the linear advection equation in one 
space dimension. With the multilevel uniform scheme, fine-grid and coarse-grid corrections 
are computed at all points of the fine grid. Then, a nonlinear filter is applied to remove the 
coarse-grid corrections at fine-grid points not contained in the coarse grid. The filtering 
produces additional errors in the form of a carrier wave with a frequency depending on 
the fine-mesh spacing. This approach does not allow for the coupling (aliasing) effects due 
to the restriction operator (fine to coarse grid transfer operator) in the multigrid method. 
However, it does offer the advantages of simplicity and application to more than two-level 
schemes. Thus, it allows the rapid comparison of multigrid algorithms. If a multigrid 
method is unstable or inefficient according to Fourier analysis of the multilevel uniform 
scheme, then it is probably not a reasonable scheme. 

In ref. 6 we consider the scalar two-dimensional advection equation and perform a 
Fourier analysis of the multilevel uniform scheme for different multigrid strategies. The 
effects of mesh-cell aspect ratio are included in the analysis. For details of the analysis, 
see ref. 6. Here, as in ref. 6, a five-stage scheme with three weighted evaluations of the 
numerical dissipation is used for a solver. The explicit stability limit of this scheme is 
extended with variable-coefficient implicit residual smoothing, which results in a Courant- 
Friedrichs-Lewy (CFL) number of 5. A two-level analysis is applied to both full coarsening 
and semicoarsening strategies. Figure 2 presents contours of the amplification factor g as 
a function of Fourier phase angles for the full coarsening and sequential semicoarsening 
strategies when the mesh-cell aspect ratio (A-R) was set to 10. Even with this AR, one 
ran clearly see the improved damping (reduced g ) in the direction of the long side of the 
cell with sequential semicoarsening. 

SPATIAL DISCRETIZATION 

A finite-volume approach, where the flow quantities are stored at the cell vertices, is 
used for the spatial discretization of the Navier-Stokes equations. For the convective flux 
calculation, an auxiliary grid is used, which is defined by connecting the cell centers of the 
original cells (see Figure 3). The inviscid numerical flux is separated into the sum of an 
averaged term that corresponds to central differencing and a dissipative term that adapts 
the discretization stencil in accordance with local wave propagation. The dissipative flux 
function is based on. the second-order-accurate upwind scheme of Yee and Harten (ref. 
7). In the case of viscous flows the entropy correction for this scheme must be carefully 
designed, as discussed in ref. 6. The physical viscous fluxes are approximated by central 
differences with a local transformation from Cartesian to curvilinear coordinates (ref. 2). 


MULTISTAGE SCHEME FOR THE FINE AND COARSE MESHES 


We have observed the need to pair spatial discretization and particular time-stepping 
schemes for the solution of the Navier-Stokes equation. The most robust choice of spatial 
discretizations found to this point is to use a second-order upwind scheme on the fine 
meshes and to set the limiter to zero everywhere on the coarse meshes. An alternative 
taken in refs. 8 and 9 is to use scalar second-difference dissipation terms on the coarse 
meshes. This approach turned out to be less robust because the second differences are 
less diffusive with respect to the acoustic modes; also, the central-difference scheme allows 
waves to travel upstream in supersonic flow. As indicated previously, a five-stage explicit 
scheme with three evaluations of dissipation is used for time advancement. Disturbances 
are most effectively expelled out of the computational domain by using local time stepping 
and implicit residual smoothing (refs. 8 and 10). The smoothing of the residuals allows the 
CFL number of the explicit scheme to be as high as 5.75, which extends the stability limit 
(CFL*) by a factor of 2.5. The time step is proportional to the ratio of the cell volume 
to the sum of the spectral radii of the inviscid flux Jacobians in the different coordinate 
directions. 

To stabilize the schemes in regions where the viscous stability limit is more restrictive 
than the inviscid limit, the coefficients of the implicit residual smoothing operator are 
locally increased, as outlined in refs. 8 and 9. At strong shocks, however, high Courant 
numbers are not appropriate. Consequently, an adaptive time step is employed. By using 
the nondimensional second difference of the pressure as a switch, the value of CFL is locally 
reduced to approximately 2 at the shock. 

MULTIGRID SCHEMES 

For the numerical solution of the Navier-Stokes equations, the two-level strategies 
presented in Figure 1 are extended to multilevel schemes, as displayed in Figure 4. The 
only differences between the two-level schemes and the multilevel schemes occur in the 
restriction process. Whenever two “down” arrows meet at a coarse mesh, averaging is 
used to obtain the restricted variable. The multilevel arrangement of the coarse meshes, 
shown in Figure 4(b), was first given by Mulder (ref. 11), who used semicoarsening to solve 
the flow alignment problem. Suitable coordinate meshes for thin boundary layers exhibit 
mostly cells with high aspect ratios in the surface-aligned direction. In this paper, other 
variants of semicoarsening, which are computationally cheaper than the semicoarsening 
schemes shown in Figure 4, are also considered for these situations. 

One may notice that the central restriction and prolongation operators discussed pre- 
viously allow for upstream propagation of disturbances in supersonic flow. Furthermore, 
the corrections given by the standard multigrid scheme near strong shocks lead to diver- 
gence of the calculation, especially when free-stream initial conditions are used. Therefore, 
the restriction operator is damped by using 

Ri,j = max (1 - (1) 

where Rij is the standard restriction operator and e itJ ^ is a switch to detect strong 
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shocks, and 


( 2 ) 


ei,/ n) = fc (n) max (vi,v i+ i,Vi-i,Vj, v j+ i,i/j-i), 


Pi— l,j 2pi ( j + Pi+l,j 

Pi-i, i + 2pi,j + Pi+ij 




Pi,j — 1 ~~ 2pt,j + Pi,j + 1 
Pi,i-i + 2pij + Pi.j+i 


( 3 ) 


where p denotes pressure. The damping coefficient is given a value of approximately 
1 in the start-up phase of the multigrid process and is decreased to a value of about 0.4 at 
later cycle numbers to allow for good asymptotic convergence rates. Such a local damping 
with a fc (n) that does not vanish is in line with the restriction damping of Karen and 
Hemker (ref. 12), who based their damping coefficients on a more physical analysis. 

A fixed V-type cycle with time stepping only on the way down is used to execute the 
multigrid strategies described above. The robustness of the overall scheme is improved 
by smoothing the resultant coarse-mesh corrections before they are passed to the finest 
mesh. The smoothing reduces the high-frequency oscillations introduced by the linear 
interpolation of the coarse-mesh corrections. The implicit residual smoothing procedure 
with constant coefficients of around 0.1 is used for this smoothing. Also, the application 
of full multigrid (FMG) provides a well-conditioned starting solution for the finest mesh 
that is considered. 


NUMERICAL RESULTS 

Two different hypersonic flow cases are used to assess the capabilities of the multigrid 
schemes. These are laminar Mach 10 ( M = 10) flow over a compression ramp and turbulent 
flow over a slender forebody at high Reynolds numbers. Table 1 gives a summary of the 
geometries and the flow parameters of the test cases. In this table, T in j is the dimensional 
free-stream temperature, and T w is the specified wall temperature. Also, the finest grid 
used for each flow computation is characterized by the streamwise and normal leading-edge 
spacings As/ e , An; e , with the normal spacing An te at the end of the geometry. 

The flow over the compression ramp is identical to case 3.2 of the Workshop on 
Hypersonic Flows for Reentry Problems, Part II, held in Antibes, France, in 1991. This 
allows comparisons with the performance of other computational methods published in 
ref. 13. Figure 5 displays the coordinate mesh generated for this test case. The low 
Reynolds number allows for a mesh with moderate aspect ratios between 5 and 50 near 
the wall. The 129 x 81 mesh is successively coarsened down to 9 x 6, which yields 9 grid 
levels with semicoarsening and 5 levels with full coarsening. The semicoarsening strategy 
is expected to eliminate most of the stiffness associated with aspect ratio. The converged 
flow solution is shown in Figure 6 for the 129 x 81 and 65 x 41 grids. The computed extent 
of separation in the corner is somewhat smaller for the coarse mesh than for the fine mesh. 
The fine-mesh results agree well with grid-converged computations published in ref. 14. 

In the next figures, we investigate the performance of the different multigrid schemes. 
For this purpose, computations were started from a solution that was converged to about 
plotting accuracy. Results from the different schemes of Figure 4 are compared in Figure 
7. The numbers indicate the final convergence rate r of the schemes and the rate of data 
processing ( RDP ) on a CRAY-YMP to advance one grid point by one multigrid cycle. 


The sequential semicoarsening scheme (Figure 4(c)) gives by far the best convergence rate. 
For this scheme, the effect of the modifications in the multigrid strategies of Figure 4 is 
investigated in Figure 8. The meshes obtained by full coarsening and by semicoarsening in 
the direction normal to the wall are both important in achieving good convergence rates. 
From Figures 7 and 8, we conclude that semicoarsening with a selected number of coarse 
meshes is most effective for this flow problem. Semicoarsening is about 2.4 times faster 
than full coarsening, which does a surprisingly good job because of its low work count. In 
particular, the multigrid scheme with sequential semicoarsening converges (i.e., the residual 
was reduced 3 orders of magnitude) in roughly 33 sec (CPU time) on a Cray-YMP, which 
is 10 times faster than the single-mesh scheme. 

The flow over a slender forebody is chosen to represent a generic configuration that 
corresponds to a high-speed civil transport aircraft or an air-breathing space transportation 
system with low wave drag. The high Reynolds numbers yield thin boundary layers, which 
can only be resolved with highly clustered coordinate meshes and large- aspect-ratio cells. 
The mesh used for the present investigations is displayed in Figure 9. The cells near the 
wall have aspect ratios up to 25000. The flow computations were done with fixed transition 
at 2 percent chord and with the assumption of an adiabatic wall. Figure 10 shows the Mach 
contours for the mesh of 256 x 96 cells, and Figures 11 and 12 show the solution obtained on 
three successively refined meshes. Both the distributions of the skin friction and the wall 
temperature are accurately computed, even with only 25 points in the normal direction. 

Next we examine the convergence behavior of the multigrid schemes. The fine mesh 
with 257 x 97 points allows 11 grid levels to be used with semicoarsening. The full diamond- 
shaped tree of coarse meshes cannot be run because the time-stepping scheme is not well 
suited to handle the extreme aspect ratios that occur on the coarse meshes. With the 
proper half of the diamond, which includes the meshes with relatively low aspect ratios, 
the numerical solution converges. Figure 13 displays a comparison of the different multigrid 
strategies. The computations are started from a preconverged solution. Again, the scheme 
with sequential semicoarsening converges best. The differences between the multigrid 
schemes for this case, which has cells with very high aspect ratios, are larger than for the 
ramp flow. The final convergence rate of the scheme with sequential semicoarsening is 15 
times better than the rate with full coarsening. A comparison of the performance for the 
complete FMG process is given in Figure 14. The sequential semicoarsening scheme takes 
194 cycles and 570 sec to reduce the averaged residuals to 10~ 2 on the fine mesh. The 
scheme with full coarsening takes 1024 cycles and 1430 sec, and the single mesh code takes 
7762 time steps and 6190 sec to achieve the same convergence level. Note that residuals of 
10 -2 correspond to a solution that is converged within plotting accuracy. If we compared 
computer times to reach lower levels of residuals instead, then the results would have been 
even better for the multigrid scheme with semicoarsening. 

CONCLUSIONS 

New multigrid schemes for hypersonic flow computations have been investigated. The 
basic solution algorithm employs upwind discretization and explicit multistage time step- 
ping. Various multigrid schemes with semicoarsening are introduced to overcome the 


486 


stiffness that results from the high-aspect-ratio mesh cells used to resolve viscous flows. 
The basic components of the algorithm are examined with a Fourier stability analysis ap- 
plied to the two-dimensional advection equation. Both the results of the Fourier analysis 
and the computations of high Reynolds number flows suggest that the semicoarsening ap- 
proach is effective. The convergence rates shown for hypersonic viscous flows are similar 
or even better than those previously published for the transonic regime in refs. 1 and 2. 
Further work is required to make the computational scheme less expensive. This need for 
more research is particularly true for the coarse meshes used within the semicoarsening 
approach, which make up the major portion of the overall work count of the scheme. 
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(a) Two levels, full coarsening, AR = 10 (CFL = 5.0, CFL* = 2.4). 
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(b) Two levels, sequential semicoarsening, weights = 1.0, AR = 10 (CFL = 5.0, CFL = 2.4). 


Figure 2. Contour plots of amplification factor for 5-stage Runge-Kutta 
scheme with first-order upwind approximation and 3 evaluations 
of dissipation (coefficients : 0.2742, 0.2067, 0.5020, 0.5142, 1.0). 
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Figure 3- Control volume for nodal-point scheme. 



(a) Full coarsening . 


(b) Semicoarsening with simple averaging. 




(c) Sequential semicoarsening- 


(d) Semicoarsening with selective averaging- 


Figure 4. Multilevel schemes. 


| ( «t mi itrnmaa 


> ■ ..hiimiw 




Figure 5. Coordinate mesh for ramp-flow problem with 128x80 cells 






Figure 7 . Influence of multigrid strategies on convergence for ramp-flow problem . 



Figure 8 . Influence of selected coarse meshes on convergence for ramp-flow problem 
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Figure 13 . Influence of multigrid strategy on convergence for forebody, mesh 256x96 . 



Figure 14 . Convergence histories for single-mesh time stepping and multigrid with sequen- 
tial semicoarsening. 
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